function mjd = frac2MJD(fd)
%frac2MJD - fractional year to modified julian date conversion
%
%       MJD = frac2MJD(FD)
%
% If FD is a vector of dates expressed in standard decimal part of a year 
% (e.g., 2-Jan-2005 is 2005+1/365 = 2005.00274), this function provides 
% the corresponding modified julian date vector MJD.
% The output date MJD is an int32 vector. For example, in the case of 
% 2005.00274, it is 53341.  
% The leap year case is automatically taken into account. 
%
% See also date2gpsw, gpsw2serial, MJD2frac, MJD2serial, serial2frac, 
% frac2serial, serial2YMD, leapy.

% G. Teza, 2021

fd = fd(:);             % to have a column vector

ud = fix(fd);
ed = fd-ud;             % Decimal part of the year

nr = size(fd,1);        % #rows   

Ily = leapy(ud);        % is 1 (logical) for leapy years
moltfact = 365*ones(nr,1);
moltfact(Ily) = 366;

nj  = round(ed.*moltfact)+1;         % "+1" is necessary to avoid the "0 day case" 

dtj = datetime([ud ones(nr,1) nj]); % datetime object

mjd = juliandate(dtj,'modifiedjuliandate');

mjd = uint32(mjd);